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Abstract 

We performed event-driven simulations of a two-dimensional granular gas between two vibrating 
walls and directly measured the local heat flux and energy dissipation rate in the stationary state. 
Describing the local heat flux as a function of the coordinate x in the direction perpendicular to 
the driving walls, we use a generalization of Fourier's law, q(x) = k\JT{x) + n\7p{x), to relate 
the local heat flux to the local gradients of the temperature and density. This ansatz accounts 
for the fact that density gradients also generate heat flux, not only temperature gradients. The 
transport coefficients k and /x are assumed to be independent of x, and we check the validity of this 
assumption in the simulations. Both k and /i are determined for different system parameters, in 
particular, for a wide range of coefficients of restitution. We also compare our numerical results to 
existing hydrodynamic theories. Agreement is found for n for very small inelasticities only. Beyond 
this region, k and ^ exhibit a striking non-monotonic behavior. 
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Driven granular gases have attracted much attention in recent years jH 0, IE B S| • This 
is partly because in these systems, energy loss in inelastic collisions is eventually balanced 
by energy input from the driving so that a stationary state can be attained. In physically 
reahstic models, the driving usually acts at the boundaries of the system, e.g. in terms of 
shearing forces or vibrating walls. Thus, maintaining a stationary state requires a subtle 
and well-balanced mechanism to transfer energy from the system's boundaries to its interior. 
The local energy (or heat) flux and the local energy-dissipation rate are therefore at the heart 
of every hydrodynamic description 0, S IE IE 1^ of th^ stationary state of driven granular 
gases. 

Fourier's law states for elastic systems that the heat current is proportional to the local 
temperature gradient, the proportionality constant being the thermal conductivity k, jlol |. 
For inelastic systems there is an additional contribution to the heat current from density 
gradients The corresponding transport coefficient /x has no analog in elastic systems. 
Theoretical approaches start from the Boltzmann-Enskog equation to account for these 
effects. Jenkins and Richman [Gj] have used Grad's moment expansion to compute the 
heat flux for small inelasticity, whereas Dufty et al. 0] have pushed kinetic models using a 
stosszahl ansatz to simplify the collision operator. MD simulations have been performed by 
Soto et al. to study /x for a granular gas on a vibrated plane in the dilute limit. 

In this paper we study the local heat flux and local energy-loss rate, as well as the 
transport coefficients k and /i for a driven granular gas in 2 dimensions. To do so, we 
perform event-driven simulations for the dynamics of N identical inelastic smooth hard 
disks of diameter a and mass m which are confined to a rectangular box with edges of 
length and Ly. Periodic boundary conditions are imposed in the ?/-direction, and the 
gas is driven through the walls perpendicular to the x-direction, see Fig. Q for a typical 
snapshot. The left and right wall vibrate in an idealized saw-tooth manner, characterized 
by the driving velocity f drive > 0: Upon a collision of a particle with the left/right wall at 
X = its x-component of the velocity changes according to —2vx ± f drive, see also 

T^ j. Inelastic inter-particle collisions are modeled using a constant coefficient of normal 



restitution a G ]0, 1[ according to n ■ = —a n ■ v^. Here, n denotes the unit vector of 
the particles' relative center-of-mass positions, and Vi2, resp. are the pre-, resp. post- 
coUisional relative center-of-mass velocities. 

A simple estimate in yields for the spatially averaged granular temperature Tq = 
N'^ X]t=i /2 in the stationary state 

2^° / 2 ^ f-) V' (l + Vl + (vr/2)%^)'. (1) 

"^^drive/^' ^ ^ 

Here e := 1 — a"^ and ip := a/2xA are dimensionless parameters. The latter involves the pair 
correlation at contact x of the corresponding elastic system and the line density A := N/ Ly. 

In the following discussion of heat flux and dissipation we will, inter alia, be interested 
in a scaling limit e | towards an elastic system. In order to prevent the system from 
heating up indeflnitely when switching off dissipation, we also need to scale the driving 
velocity fdrive =: ^Vq, where Vq is flxed. Hence, the driving vanishes in the elastic scaling 
limit, and the spatially averaged temperature reaches a flnite value, approximately given by 
To {2/TTf{xX)~'^mvQ according to Eq. ((T]). 

Locally, the translational energy changes due to collisions as well as due to free streaming 
of the particles in between collisions. The latter gives rise to a kinetic contribution to the 
heat current, qr"^™ , whereas the collisions are responsible for the energy loss ( as well as for 
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FIG. 1: Model of N disks, driven in the x-direction with periodic boundary conditions in the 
y-direction. 



the collisional contribution to the heat current q™*. For vanishing macroscopic velocity the 
energy balance equation is usually formulated as 0| 

p(r, t)^T{r, t) = -V ■ q{r, t) + C(r, t) (2) 

for the hydrodynamic fields of density p, temperature T, total heat current q = q^^"^ + q'"* 
and local energy dissipation rate (. 

To compare our simulations to the hydrodynamic theory, we need to introduce a coarse 
graining function $(r) which is nonzero in a small area centered at r only. We require 
of course J dr^{r) = 1. In the absence of a velocity field, the coarse grained kinetic heat 
current is defined by n 2 

q^-{r,t)=J2^v,^r-r,) . (3) 
1=1 

The coarse grained density p{r) and temperature T(r) are defined analogously. 

The change of energy during a binary collision in a small time interval At can be decom- 
posed into a source term and a divergence of a flux. The coarse grained energy dissipation 
rate C{fyt) can be calculated analogous to [l^ and is given by 

in terms of the change of energy AEi\j of particle i due to a collision with particle j in 
the time interval [t,t + At]. The prime at the summation sign restricts i and j to those 
particles colliding in At. The energy dissipation rate trivially vanishes in the elastic case, 
when AEi\j = —AEj\i. Similarly, the collisional contribution to the heat current is given by 
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X / $(r - ri + svij) ds, (5) 
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FIG. 2: The heat flux in the x-direction for a system of = 256 particles in a box of size 
Lx = 20 and Ly = 25, corresponding to a total area fraction of (f)Q = 0.4. The dotted line shows 
the kinetic part g^™, the dashed line is the collisional part gj?*, and the full line represents the total 
heat flux q^- The inset displays the local energy-loss rate Q. 

where r^j = — Vj. In our simulations we have never found any significant variations of 
the long time averages of the hydrodynamic fields in the y-direction parallel to the driving 
walls. Hence we coarse grain our system by subdividing the box into stripes of width Ax. 

In the following we report numerical results on the stationary state only. The reader who 
is interested in more details of the simulations, such as initialization or relaxation towards 
the stationary state, is referred to jl^]. For dimensional reasons, the driving velocity Vq sets 
the time and energy scale and is chosen to be 1. Likewise, the particle mass and diameter 
set the mass and length scales and are set to 1, too. 

We first present data for the heat-flux profile in Fig. El for a system of global area fraction 
00 := NTr/{4LxLy) = 0.4 and coefficient of normal restitution a = 0.9. More precisely, the 
figure displays the x-component of the heat flux and its kinetic and collisional part as a 
function of the normalized position x/L^. The inset shows the local energy-loss rate and 
will be discussed below. The heat flux is antisymmetric about the middle of the system as 
expected. We clearly see that its collisional part ©, which is represented by the dashed 
line cannot be neglected. It is of the same order as the kinetic contribution — depicted by 
the dotted line and corresponding to Eq. (jS)). For lower density systems the collisional 
contribution becomes less important but still it is not negligible for global area fractions 
even as low as 0o = 0.1. 

Now we turn to the dependence of the heat flux on the coefficient of restitution a in the 
elastic scaling limit. We recall that the driving strength has been adjusted such that the 
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FIG. 3: Rescaled heat flux (1 — a^) ^(1 + a) '^qx for a wide range of coefficients of restitution 
0.6 < a < 0.995 for otlierwise fixed systems (L^ = 20, N/Ly = 10.24 [0o = 0.4]). 

granular gas attains a finite temperature as e —>■ 0. The heat flux is proportional to the 
temperature per time and hence expected to scale like Tq f drive ~ cf. (0). This argument 
is checked by plotting q^/e for different coefficients of restitution. The collapse of different 
data sets is excellent for almost elastic systems (e <^ 1). With a slight modification of the 
scaling according to qx/[e{l + a)^], the data collapse also works approximately for the full 
range 0.6 < a < 0.995 and is shown in Fig. El 

In Fig. m we show the heat flux for a fixed coefficient of restitution a = 0.99 and various 
system sizes L^. It has been shown that the scale for temperature inhomogeneities is 
set by Lx, so that we expect the heat flux to scale like L~^, if plotted versus x/Lx- In Fig. HI 
we plot qxLx for various system sizes Lx and get a decent data collapse for dilute systems. 
For higher densities this scaling captures at least the correct order of magnitude. 

It is a characteristic feature of hydrodynamics of inelastically colliding particles that a 
heat current can be generated not only by a nonuniform temperature but also by density 
inhomogeneities. If the spatial variations of temperature or density are restricted to long 
wavelengths, one would expect a gradient expansion to hold. The simplest constitutive 
equation for the heat flux is thus a straightforward generalization of Fourier's law, as dis- 
cussed in the literature, cf. |lfl|]. Chapman-Enskog expansions of both the Boltzmann and 
the Boltzmann-Enskog equations predict that the heat flux of an inelastic system is given 



where k is the heat conductivity and is a new transport coefficient that has no analog for 
elastic systems 0| • 



by 



— T(x) + — p(x) 



(6) 
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FIG. 4: Rescaled heat flux L^Qx for a wide range of box edges (area fraction 0.001 < (f)o < 0.4 
at fixed line density N/Ly = 10.24 and fixed coefficient of restitution a = 0.99. This graph looks 
almost the same for all a > 0.9. 



To estimate the transport coefficients from our data, we assume that both, k and yU, do 
not depend on the position x, i.e. we assume them to be constant throughout the box. It is 
then straightforward to extract them from a fit of our data to the above expression In 
Fig. we show both transport coefficients k and ix/e. Both are non-monotonic in a. In the 
elastic limit we find that k tends to a non-zero constant, while n oc e. The fit of the heat 
fiux as computed from Eq. (jU)) with constant k, and fi is very good for all investigated a, e.g. 
for a = 0.9 it cannot be distinguished from the data shown in Fig. [21 

Both transport coefficients have been computed within kinetic theory for small inelastic- 
ities by Jenkins and Richman |^, who find for the thermal conductivity 



2 + 30x^2(4r - 3) 
rx(17 — 15r) 



(2 + 30xr) + 



(7) 



where r := (1 + a)/2. We use the Henderson approximation for the pair correlation at 
contact X which is a function of the area fraction only. Inserting the local temperature T{x) 
and local area fraction from our simulations into ((7j), we obtain a spatially dependent 
transport coefficient k\oc{x) and, from a corresponding equation in @], fj^iodx). The resulting 
heat fiux from (jH)) is shown in the inset of Fig. E] and found to agree very well with the 
simulations for a > 0.99. The agreement is reasonable even up to a ~ 0.96. For larger 
inelasticities the theoretical curves capture the correct order of magnitude, but overestimate 
the curvature of qx{x). 

Alternatively we can use the global or mean temperature Tq and area fraction 0o to 
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FIG. 5: Transport coefficients k and fi as functions of the coefficient of restitution a for systems 
of size Lx = 20 and line density N/Ly = 10.24 (0o = 0.4). 

evaluate Eq. (|7j) and the respective equation from for /i. In Fig. El we compare this 
to the thermal conductivity k, we obtained from fitting our simulations to Eq. The 
agreement is good as long as a > 0.99. As to fio (not shown), the deviations to /x are quite 
strong. Note, however, that the theory of Ref. p as well as the simulations of Ref. 
apply to the low density limit and not to 0o = 0.4. Furthermore, the difficulties in our 
fitting procedure increase considerably as a — 1 because the temperature gradient becomes 
proportional to the density gradient. Consequently, it is not possible for « — 1 to determine 
two parameters from the fit unambiguously. 

The difference between k,\oc{x) and kq is only a few percent for a > 0.99. The difference 
in heat flux is even less, because the strongest inhomogeneities in the transport coefficients 
occur in the middle of the sample where the gradients of temperature and density vanish. 
The heat flux, computed using the constant transport coefficients hq and fiQ, is also shown 
in the inset of Fig. 0] for comparison. The fiuctuations of the transport coefficients with x 
increase with increasing inelasticity, e.g. for a = 0.9 we find kioc(0)/k ~ 1.5. 

To estimate the degree of inhomogeneity we have divided the box into an inner and an 
outer part and fitted the heat current using data from either half of the box only. The 
scattering of the data from the inner and outer part is shown in Fig. Eland provides a rough 
measure for the effects of inhomogeneous transport. The results for the outer half are almost 
identical to the ones for the full system. This is also true for the inner part except for the 
weakly inelastic systems for which the absolute value of the heat flux in the middle of the 
sample is so small that statistical fluctuations dominate. 

The local energy loss, as deflned in Eq. (0]), is shown in the inset of Fig. Again, 
for a > 0.99 the agreement with the predicted C = —IQecf^xiT / 'kY^'^ from |(|| is very 
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good (not shown). For very dilute quasi-elastic systems (0o ^ 0.1, a > 0.95) the dissi- 
pation is highest in the middle of the sample. For denser and/or more inelastic systems 
we find the absolute value of C, (the dissipation) to be greatest in an intermediate region 
X ~ ±0.3 L3.. Even though the density is largest in the middle of the sample, energy dissipa- 
tion is not very pronounced there because the mean kinetic energy is already comparatively 
small. Integrated over the whole box, the local energy loss has to fulfill the conservation law 
qx{—Lx/2) — qx{+Lx/2) = /^2''/2 Ci^)dx. This, of course, is confirmed in the simulations, 
but requires careful measurements of the heat fiux at the boundaries. 

In this paper we have discussed the heat current in a granular gas driven by vibrating 
walls. We have measured the collisional contribution to the heat current and have shown 
that it cannot be neglected, not even in low density systems. We have extracted transport 
coefficients, the thermal conductivity k, as well as the transport coefficient fi, which accounts 
for a heat fiux due to density inhomogeneities. Both transport coefficients have been deter- 
mined for a wide range of inelasticities 0.6 < a < 0.995 and area fractions 0.01 < (po < 0.4. 
We have studied in detail the elastic limit which is reached by scaling the driving velocity 
with e, so that the temperature tends to a finite value as a — 1. Comparing our data to 
theoretical work in two dimensions |6|] we found good agreement for a > 0.99. For stronger 
inelasticities non-monotonic behavior is observed: The thermal conductivity k first decreases 
with increasing e, goes through a minimum around a ~ 0.96 and then increases again. The 
transport coefficient fi/e first increases with e, goes through a maximum at a ~ 0.96 and 
then decreases again. Furthermore a rough estimate of the fiuctuations of the transport 
coefficients with spatial position in the box has been given. These effects are expected to be 
strong for moderately inelastic systems which are characterized by nonuniform density and 
temperature and need further investigation. 

We thank Hans Vollmayr for pointing out to us the importance of the elastic limit and 
Isaac Goldhirsch for suggesting to measure the heat flux. We acknowledge flnancial support 
by the DFG through SFB 602, as well as Grant Nos. Zi 209/6-1 and Mu 1056/2-1. 
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